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ABSTRACT 

£S) . We derive the properties of dusty tori in Active Galactic Nuclei (AGN) from the 

' comparison of observed Spectral Energy Distributions (SEDs) of SDSS quasars and 

a precomputed grid of torus models. The observed SEDs comprise SDSS photometry, 
■ 2MASS J, H, and K data, whenever available and mid-Infrared (MIR) data from the 

S pitzer Wide- a rea In fraRed Extragalactic (SWIRE) Survey. The adopted model is that 
of iFritz et al.l (|2006t ) . The fit is performed by standard \ 2 minimisation, the model 
£N| ' however can be multi-component comprising a stellar and a starburst components, 

' whenever necessary. Models with low equatorial optical depth, T9.7, were allowed as 

well as "traditional" models with T9.7 ^ 1.0, corresponding to Ay ^ 22 and the results 
were compared. Fits using high optical depth tori models only produced dust more 
compactly distributed than in the configuration where all T9.7 models were permitted. 
Tori with decreasing dust density with the distance from the centre were favoured 
while there was no clear preference for models with or without angular variation of 
the dust density. The computed outer radii of the tori are of some tens of parsecs large 
but can reach, in a few few hundreds of parsecs. The mass of dust, Moust, 

and infrared luminosity, Lir, integrated in the wavelength range between 1 and 1000 
Jim, do not show significant variations with redshift, once the observational biases are 
taken into account. Objects with 70 um detections, representing 25% of the sample, 
are studied separately and the starburst contribution (whenever present) to the IR 
luminosity can reach, in the most extreme but very few cases, 80%. 
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1 INTRODUCTION 

Efforts in understanding the physics at work in active galac- 
tic nuclei (AGN) started four decades ago. The origin of the 
infrared (IR) continuum of AGN was initially a matter of 
controversy, as it could be non-thermal but could equally 
be due to ther mal emission fro m dust grains. It was long 
ago suggested |Rees et al.lll96sl ) that IR emission radiation 
from Seyfert galaxies in the 2.2 - 22 um wavelength range 
was produced by dust grains heated by ultraviolet (UV) and 
optical emission from the nucleus. Work carried out later 



on suggests the IR emission to be the reprocessed emis- 
sion of the UV /optical radiation from the accretion disk 
by the particles co mposing the torus, namely silicate and 
graphite grains (e.g.lPier fc Kr olik 1992 ; iGranato fc Danese 



1994 : lEfstathiou fc Rowan-Robinsonl 1 19951 : iNenkova et al 
20021 ). Various configurations of the dust distribu- 
tion geome try and compositions have been since sug- 



geste d (e.g. Pier fc Kroliklll992l: Ivan Bemmel fc Dullemondl 
20031: iDullemond fc van Bemmell 120051 : IFritz et all 1200a 

Elitzur fc Shlosman|[20ul K 

Recent observations (|jaffe et al indicate that this 
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structure might be smaller than originally thought, but re- 
solved in at least some nearby systems. Interferometric mea- 
surements at 8 - 13.5 pm centred on NGC 1068, a prototype 
Seyfert 2 galaxy, revealed a structure wit h an extent of ~ 3 .4 
pc that can be identified with the torus (|jaffe et al.ll2004 ). 

Dust remains essential to understanding the Seyfert 1/2 
dichotomy. The dust sublimation radius defines the outer 
boundary of the broad- line region (BLR). Recent results 
from reverberation m apping place the IR-em itting medium 
just beyond the BLR l|Suganuma et alj|2006h . The IR emit- 
ter presumably corresponds to the torus, that would make a 
type-1 object look like a Seyfert 2 nucleus when the line of 
sight to the BLR is obscured by the dusty medium. However, 
if the obscuring torus is "clumpy" ra t her th an homogeneous, 
as s uggested by Krolik fc Begelmanl l| 19881 ) and modeled by 
e.g. iNenkova et al. 

I (|2002T ). obscuration and classification of 
AGN becomes probabilistic. 

In this paper we address the question of the distri- 
bution of properties of dust tori around AGN, studying 
SDSS quasars in the fields covered by the Sp itzer Wide-area 
InfraRed Extragalactic Survey (SWIREo Lonsdal e et al.l 
2003). We assume that tori are smooth (as opposed to 
"clumpy") distributions of dust, composed of graphite and 
silicate grains and will limit ourselves to the use of pho- 
tometry only, spanning a large wavelength range, from ~ 
0.35 to 160 pm. Our aim is to classify the properties of 
the objects under study fitting their observed Spectral En- 
ergy Distributions (SEDs), built by photometric points col- 
lected from various surveys, to various model components, 
including stellar templates, AGN tori and starburst emis- 
sion, whenever applicable. The quasar sample consists of all 
spectroscopically confirmed type-1 quasars in the common 
regions between SWIRE and S PSS Data Release 4 (DR4; 
lAdelman-McCarthv et~al1l2006h . The purpose of this work 
is to constrain the model parameters and to quantify the IR 
properties of bright quasars, within the limitations imposed 
by the model assumptions, the resolution provided by the 
broad band photometry, and the degeneracies resulting from 
the fitting procedure. 

The paper is structured as follows. Section [5] describes 
in detail the samples. A brief description of the torus models 
and the stellar and starburst components is given in Section 
[3] Section [3] describes the SED fitting mechanism and the 
physical parameters thus derived. The results are presented 
in Section [5j with a subsection dedicated to the problems 
of degeneracy and aliasing and their implications. Finally, 
Section [6] presents a discussion of the results of this work. 



2 THE SAMPLE 

For the purposes of this study consistent IR data are es- 
sential, covering as large a wavelength range as possible. 
SWIRE provided the astronomical community with un- 
precedented quality MIR photometric data for over 2 mil- 
lion objects, including hundreds of thousands of AGN of 
all types. We therefore selected the samples from regions 
with overlapping SDSS DR4 and SWIRE IRAC or MIPS 
24 pm data; SWIRE ELAIS Nl and N2 and the Lockman 



field. The sample comprises 278 spectroscopically confirmed 
SDSS quasars within these fields with redshifts spanning 
0.06 < z < 5.2. 

The SED of each object was constructed using the SDSS 
photometry (u, g,r,i, z), 2MASS photometry ( J,H,K) and 
IR photometry from SWIRE, i.e. IRAC 1-4 channel fluxes, 
MIPS 24, 70 and 160 urn fluxes, whe never available. 2MASS 
photometry from the 2MASS x6 (jBeichman et al.l l2003h 
was used, for objects in the Lockman hole. The SWIRE 
catalogues we use throughout this work were processed by 
the SWI RE collaboration . Deta il s about the d ata c an be 
found in lLonsdale et~ai1 (|2004h . ISurace et ail (|2004T ) and 
IShupe et al.l ( 2007 ) . The selection of "reliable" sources in 70 
and 160 pm is done as follows. All sources with detections 
above 5<r are considered reliable. The 70 and 160 pm flux 
limits adopted correspond to 90% completeness and gen- 
erally coincide with the 5a noise of the images (measured 
from the sky rms). They slightly vary with the field, and for 
SWIRE ELAIS Nl and N2 and the Lockman fields they are 
of 17, 17.5 and 18 mjy for the 70 micron, and 104, 124 and 
108 mjy for the 160 pm . Also sources with at least 3<r de- 
tections and a 24 micron counterpart brighter than 300/iJy, 
within 6" and 12", for the 70 and 160 micron sources, are 
taken to be real. There will be, therefore, cases where the 
70 and 160 pm fluxes will lie below the nominal flux limits. 
For more details, se e Vacc ari et al. (in prep). 

iRichards et al] l|2006h presented an extensive multi- 
wavelength analysis of a sample of 259 objects, most of which 
belong to our quasar sample, also. For the colour properties 
and global SEDs, we therefore refer to this work. For a de- 
tailed analysis of the properties of the quasars and composite 
SEDs in SWIRE EN1 field, see also lHatziminaoglou et al.l 
l|2005h . 

Tables [T] and [2] list the coordinates, redshifts, optical, 
near- and mid-IR photometry of the 278 quasars. Here only 
the first 20 entries are presented, the full tables are available 
as on-line material. The SWIRE fluxes missing from Table 
[2] do not indicate drop-outs but objects lying at the edges 
of the fields and escaped detection in some of the SWIRE 
bands due to slight variations in the rastering. As only a 
relatively small fraction of the sources have been detected 
in 70 pm, the 70 and 160 pm photometry will be presented 
in a separate table in Section f5. 51 



http:/ /swire. ipac.caltech.edu 



3 THE MODEL COMPONENTS 

The observed SED of a galaxy, from the UV to the far-IR 
(FIR), can in principle be built as the sum of three distinct 
components: stars, which emit most of their power in the 
optical, hot dust, heated by accretion onto a supermassive 
black hole, whose emission peaks in the MIR (from a few to 
some tens of microns) and cold dust, mainly heated by star 
formation. 



3.1 The torus component 

The main focus of this work is the dusty torus, for which we 
are seeking to determine the properties for the individual 
objects and for the sample as a whole. For the purposes of 
this work we assume a continuous dust distribution around 
the central sources (an accreting black hole), consisting of 
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Table 1. Coordinates, redshifts and SDSS ugriz photometry for the first 20 quasars, ordered by RA. 
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Table 2. 2MASS photometry, IRAC (in fijy) and MIPS 24 pm fluxes (in mjy) for the first 20 quasars, ordered by RA. 
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silicate and graphite grains, confined in a toroidal shape, 
as op posed to "clumpy" tori models (see e.g. iNenkova et al.l 
2002). For a deta iled description of the model we refer to 
Fritz et all l|2006l ). Here, however, we summarize some of 
the major characteristics of the model. 

The central source is assumed to be point-like and its 
emission isotropic. Its spectral energy distribution is denned 
by means of a composition of power laws (i.e. XL(X) oc A 1 ) 
with different values for the spect ral index i. We adopted 
va lues found in the literat ure (e.g. iGranato fe Danesdll994l 
or lSchartmann et al.ll2005l 'l. namely: i = 1.2 for 0.001 < A < 
0.03, i = for 0.03 < A < 0.125, i = -0.5 and i = -1.0 for 



0.125 < A < 10, and a Rayleigh- Jeans decline (i — —3.0) is 
assumed longward of A = 10, where A is given in pm . 

The torus geometry adopted to describe the shape and 
the spat ial distribution of dust is the so-called "flared disk" 
(see e .g. lEfstathiou fc Rowan- Robinson 1995, Mansk e"et al.l 
Il998l and Ivan Bemmel fc Dullemondll2003l ). that is a sphere 
with the polar cones removed. Its size is defined by its outer 
radius, K ou t, and the opening angle, O, of the torus itself. 
The dust components that dominate both the absorption 
and the emission of radiation are graphite and silicate. The 
location of the inner radius, R in , depends both on the subli- 
mation te mperature of the dust grains (1 500 and 1000 K, for 
graphite (|Barvainid Il987l ) and silicate (jGranato fc Danesd 
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Il994h . respectively) and on the strength of the accretion 
luminosity. We adopted the absorptio n and scattering coef- 
ficients given bv lLaor fc Draind (|l993h for dust grains of dif- 
ferent d imensions, weighted with the standard MRN distri- 
bution |Mathis et al.|[l977h . Grains dimensions range from 
0.005 to 0.25 um for graphite, and 0.025 to 0.25 um for sili- 
cate. 

The dust density within the torus is modeled in such 
a way to allow a gradient along both the radial and the 
angular coordinates: 

p(r,0) = a .Ae- 7|cos(fl)l (1) 

with (3 taking the discrete values 0.0, -0.5 and -1.0, and 7 the 
values 0.0 and 6.0. When 7 7^ 0.0 the dust distribution in no 
longer a "flared disk" but takes a shape that resembles that 
of a donut, hence the name "torus". The "zero value" a is 
defined by the value of the equatorial optical depth at 9.7 mi- 
cron, Tg.7. One of the novelties of this work is the use of low 
optical depth tori (rg.7 < 1). Even thou gh they are part of 
mode l sets present in the literature (e.g. iGranato fc Danesei 
Il994h . they have not been used to explain the IR emission 
of AGN. The implications of this approach are explained in 
Section [5] 

The global model SEDs are computed at different an- 
gles of the line-of-sight with respect to the torus equatorial 
plane, in order to account for both type-1 and type-2 objects 
emission, and include three contributions: emission from the 
AGN (partially extinguished if the torus intercepts the line- 
of-sight), thermal and scattering emission by dust in each 
volume element. 

Torus models with R ou t /Rin of 300 are a priori excluded 
from the runs, even though they belong to th e "standard" 
grid of models presented bv lFritz et al.l l|2006h . as they im- 
ply tori with physical sizes of several hundred parsecs, some- 
times even kpc. For such a ratio, an AGN of 10 46 erg/sec 
accretion luminosity and an inner radius of 1.3 pc would 
have an outer radius of ~400 pc. In the general case, only 
models with R ou t /Rin of 30 and 100 are allowed, but models 
with Rout /Rin of 300 will be revisited in Section[53]in order 
to address specific cases. 

3.2 The stellar component 

The stellar emission will play a minor role in this work, since 
we are dealing with bright, mainly high redshift quasars, for 
which the galaxy light is not significant. However, for com- 
pleteness, we add a stellar component modeled as the sum of 
spectra of Simple Stellar Population (SSP) models of differ- 
ent age, all assumed with a common (solar) metallicity. The 
set of SSP used for this analysis is b uilt with Padova evo- 
lutionary tracks (iBertelli et al.lll994 ). a Salpe ter IMF with 
masse s in the range 0.15 - 120 Mq and the Ijacobv et al.l 
(1984) library of observed stellar spectra in the optical do- 
main. The extension to the UV and IR range is obtained by 
means of Kurucz theoretical libraries. Dust emission from 
circumstellar envelo pes of AGB stars has been added by 
iBressan et al] i|l998ft . 

3.3 The starburst component 

For the cold dust component, the major contributor to the 
emission at wavelengths longer than ~ 30 um, we use two 



observational starburst templates, namely M82 as a repre- 
sentative of a "typical" starburst IR emission and Arp 220 
as representative of a very extinguished starburst. A more 
exhaustive approach would require us to provide a physical 
model of the starburst component, which however is far be- 
yond the scope of this work. Starburst templates are used 
only when there are observed datapoints longward of 24 um 
restframe, (typically only when 70 and/or 160 pi data are 
available) , and if a torus model fails to provide an acceptable 
description of the observational SED. 

This choice is, in fact, arbitrary as nothing forbids the 
presence of a starburst component fainter than the 70 um 
detection limit. The flux at 24um, though, is dominated by 
the torus: the starburst contribution at this wavelength is 
minimal, as it coincides with the presence of a deep absorp- 
tion feature in their SEDs (see Section [53}, an with no more 
datapoints redwards of 24 um it is simply impossible to con- 
strain that part of the SED. Adding a starburst component 
in order to fit all objects would only increase the degeneracy 
and uncertainties of our results. Quantities such as the IR 
luminosity (see Sections 14. Il and l5.2|l . though, might be seen 
as lower limits for the objects with no MIPS 70 and/or 160 
um detections. 



4 SED FITTING 

Given the large amount of data available, we developed a 
fully automatic fitting procedure, where the goodness of the 
fit is measured in terms of a x function: 

where Oi are the observed values, Mi the values computed 
from the model and at are the observed errors of the i-th 
photometric point. The expected values from the model are 
computed by convolving the synthetic flux with the filters' 
response curves, after an opportune normalisation and K- 
correction is applied. The dominant component in the UV 
and NIR (rest-frame) is the accretion disk. For very low 
redshift objects light from the host galaxy might also be 
present. The former is clearly distinguishable from a typical 
stellar SED, because it is in general bluer and flat over the 
entire range of overlap. Hence, if a good fit is not achieved, 
stars are removed from the final SED, and a pure AGN com- 
ponent is used at these wavelengths. Since we are dealing 
with an AGN sample (i.e. we expect an AGN component 
to be present in the observed SEDs of all the objects), the 
starburst emission, which dominates over the other compo- 
nents only in the FIR, is included only if there are observed 
70 and/or 160 um datapoints. 

Examples of fits are shown in Fig. [I] on the left, an 
object whose SED was reproduced by an AGN component 
only; in the middle a case of an AGN with starburst emis- 
sion; and in the right an object with all three components 
(torus, starburst and stars) present. 

Even though the minimum \ 2 wm define the best fit, 
the associated probabilities can not be taken at face value, 
as for a number of reasons the derived reduced y^s are over- 
estimated. First of all, in order to compute the model mag- 
nitudes, the models are convolved with the filters' transmis- 
sion curves. If the model is a very accurate representation of 
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Figure 1. Examples of fits: an object whose SED was reproduced by an AGN component only (left); a case of an AGN with starburst 
emission (middle); and an object with all three components (torus, starburst and stars) present (right). 



the real SED of an object, this convolution will lead to more 
accurate results. In our case, however, the optical/UV part 
of the model is simply a power law, as already mentioned in 
Section 13.11 while the SDSS photometry, that corresponds 
to the same part of the SED is quite sensitive to the pres- 
ence of broad emission lines or the presence of the small blue 
bumpE Also, as a general remark, the photometric errors are 
very small (typically of the order of few percent, as seen in 
Tables Q] and for both SDSS and SWIRE datapoints. And 
even though we use the errors in the catalogues to properly 
weight the fits, in many cases the computed values of the 
reduced x 2 s are very high due to the small photometric er- 
rors. In order to avoid excessively high weighting and high x 2 
values, many SED fitting codes imp ose mini mum flux errors 
(e.g. HyperZ jBolzonella et al.ll200ol ; ImpZ. iBabbedee et al.l 
2004), here however we chose not to adopt this approach. 

4.1 From SED fitting to physical parameters 

Based on the best-fit model parameters a number of other 
physical parameters can be derived: 

• the accretion luminosity, L acc ; this is, the soft X-ray, 
UV and optical luminosity coming from the accretion disc, 
which provides the main source of dust heating in the torus 

• the IR luminosity, Ltr, defined as the integral of all the 
components in the interval 1 to 1000 um 

• the relative contribution of AGN and starburst activity 
to the IR luminosity (the latter is computed only where 70 
and/or 160 um data are available), i.e. the fraction of each 
component with respect to the total IR luminosity 

• the innermost radius of the torus, Ri n , which is the dis- 
tance at which the grains reach their sublimation tempera- 
ture, averaged over all sizes of graphite grains: 

Rin~ 1.3- VWTfsM M, (3) 

where T1500 is the sublim ation temperatu re of the dust grain 
given in units of 1500 K (|Barvainislll987r ) 

2 The optical/UV continuum, however, is very well constrained 
by the available datapoints and therefore the computation of the 
accretion luminosity is not affected - see also Section |4. II 



• the optical depth (or extinction in the V band) along 
the line of sight 

• the hydrogen column density along the line of sight 
(note that in this case the Galactic dust-to-gas ratio is im- 
plicitly assumed, as a consequence of the use of the MRN 
distribution function) 

• the torus full opening angle; this parameter also defines 
the covering factor, which is the percentage of the solid an- 
gle which is covered by dust in the torus, as seen from the 
nucleus 

• the mass of dust, Mnust, within the torus, which is the 
integral of all dust grains over all volume elements 

4.2 High and low optical depth, T9.7, models 

We also address the issue of low optical depth tori, namely 
tori with equatorial T9.7 1.0, equivalent to Av ^ 22. The 
obscuring medium (torus) is usually considered to be opti- 
cally thick but there are no physical arguments against the 
existence of optically thin tori. We therefore test the possi- 
bility of quasars seen also through low optical depth tori, as 
opposed to the "traditional" picture of them seen uniquely 
on lines of sight not intercept by the torus. In order to do so, 
we run the SED fitting twice, once allowing for all optical 
depths, and once allowing only for models with (rg.7 ^ 1-0), 
and compared the results. 



5 RESULTS 

From the 278 quasars comprising the sample, 247 have 
IRAC coverage, 86 of which have additional JHK data from 
2MASS (of a total of 97 quasars of our sample detected by 
2MASS) and a total of 268 has a 24 pm datapoint. If one 
requires a good sampling of the SEDs, one might require de- 
tections in at least 8 bands: u,g,r,i,z, two out of four IRAC 
bands (due to slight differences in the rastering, objects that 
lie at the edges of the fields might escape detection in 2 of 
the 4 bands) and MIPS 24 pm. Similar considerations ap- 
ply in that the objects that were not detected in this band 
were not dropouts but simply happened to lie outside the 
covered areas. This condition is fulfilled by 237 of the 278 
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Figure 2. Minimum \ 2 distribution in the case all models are 
allowed in the fits (red plain line) and that of only high rg.7 tori 
models allowed (black dashed line). 

objects. For 70 of these there are 70 pm detections. The 41 
objects that do not have IRAC data or 24 um detections 
were included in the sample but it will be made clear when- 
ever their exclusion causes changes in the results mentioned 
below. 

Figure [2] shows the minimum \ 2 distribution in the case 
all models are allowed in the fits (red plain line) and that 
of only high rg.7 tori models allowed (black dashed line). 
For the statistical analysis of the sample and based on the 
discussion presented in Section f3] we consider the objects 
for which the best fit solutions have reduced % 2 lower than 
16 for both runs. This cut is of course arbitrary and it has 
been chosen because it corresponds to a minimum in the \ 2 
distribution corresponding to the run where all rg.7 where 
allowed (see Fig. [2} . This cut will exclude in both runs 10% 
of the objects. Note that the x 2 is narrower in the case where 
all models are allowed, but broadens when constraints on the 
model grid are imposed. Having chosen a more conservative 
cut in the values of reduced x 2 , e -g- 10, would have excluded 
another ~4% of the fits and would not further affect the 
results. 

From the 250 of the 278 with reduced x 2 < 16 , 190 were 
better matched with an AGN component alone, 46 required 
an additional starburst component, five more were sets of 
AGN, starburst and stellar component, and another nine 
were just composed of stars and an AGN. 

Tori with j3 equal to -1.0, were slightly favoured, with 
40% of the observed SEDs being reproduced by such models, 
with /3=-0.5 and 0.0 providing better fits to 32% and 28% 
of the objects, respectively. Models with 7 equal to 6.0 (see 
equation [TJ provided better fits to 46% of the sample while 
there was a slight tendency towards tori with R ou t /Rin of 100 
(59% of the objects). However, when high Tg.7 were forced 
to the fit, the tendencies were reversed, with 56% of the ob- 
jects favouring models with 7 equal to 6.0 and R out /Ri n of 30 
and 50% of the objects better matching (3 = —1.0. All these 
changes tend to produce dust more compactly distributed 
than in the configuration where all Tg.7 models were per- 
mitted and this change occurred most probably in order to 
avoid large amounts of IR emission to be produced otherwise 
by the model. 

Fig. [3] shows the distribution of the best values for the 



Figure 3. Distribution of viewing angle (in degrees), from the 
equatorial plane. The red solid line (black dashed line) corre- 
sponds to the run with all tori models allowed (only high equa- 
torial optical depth, rg.7, models allowed). 

viewing angles, measured from the equatorial plane towards 
the pole, for objects with minimum x 2 s below the selected 
threshold. The colour coding is the same as in Fig. [2] 

When all optical depths are allowed, the SEDs of more 
objects could be reproduced with configurations allowing 
for the line of sight (LoS) to intercept a (low rg.7) torus. 
When only high rg.7 values are allowed, the histogram of 
viewing angles is depleted towards the low angles end and 
boosted towards the large angles end, as expected from our 
current view of the AGN paradigm, in which type-1 objects 
can not be seen through the torus. In both cases, however, 
for the majority of objects the viewing angle is larger than 
half the torus effective angle, and therefore the LoS does not 
intercept it. 

5.1 The covering factor 

The covering factor (CF), as already mentioned, is the per- 
centage of the solid angle blocking the light of the AGN due 
to dust absorption. For models with high equatorial optical 
depths, rg.7, it is proportional to the cosine of the effective 
opening angle of the torus. Fig. f3] shows the distribution of 
CF for the run with all models allowed and that for which 
only high values of rg.7 were used. The gaps are due to the 
discrete values given to the torus half opening angle when 
the grid of models was built, namely of 20, 40 and 60 degrees. 
For low rg.7 models and models with density decreasing to- 
wards high altitudes above the equatorial plane (7 > in 
equation [1} though, the covering factor has to be recom- 
puted to account for the fact that the torus obscuration is 
considered effective only where r(A = 0.3/i) > 1. 

The colour coding is the same as that in Fig. [2] The plot 
implies that very low optical depth tori correspond to very 
small covering factors (first bin of red solid histogram). This 
happens because in these cases the CF is recomputed to ac- 
count for the fact that the quasars could be seen through an 
optically very thin medium. We consider that the LoS inter- 
cepts the torus when the optical depth at 3000A reaches the 
value of 1. When these models are excluded, the CF assumes 
a steeper distribution, with many more objects having much 
larger dust coverage. 
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Figure 4. Distribution of the best fit values for the torus CF. 
Colour-coding is the same as for Figs. [2] and [3] The gaps are due 
to the discrete values given to the torus half opening angle when 
the grid of models was built (see text fro more details) . 
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Figure 5. L aC c/LiR as a function of the viewing angle and the 
accretion luminosity, L acc - 



5.2 IR luminosity and mass of dust 

For objects whose MIR SED can be represented by a torus 
component only, the ratio of L acc over the IR luminosity, 
Lir, would be another indicator of the obscuration. Infrared 
radiation is UV/optical radiation reprocessed by the dust 
composing the torus. For type-1 objects, for which the cen- 
tral engine is directly observed, this ratio is expected to take 
values larger than unity, as indeed is the case. Furthermore, 
Lacc /Lir does not show any dependence neither on redshift 
nor on L acc - There is, however, a clear tendency of decreas- 
ing L a cc /Lir with increasing viewing angle shown in Fig. [S] 
with a factor of 4 of difference between 10° and 90°. This ten- 
dency is due to the fact that Lir is orientation-independent 
while the measurements of L acc are not, translated into a 
constant Lir and a decreasing L acc with increasing view- 
ing angle. From the results of the best fit models one notes 



Figure 6. Mass of dust, Mousti as a function of redshift, z, and 
accretion luminosity, L acc . 



that objects with best viewing angles of less than 30°, are 
also assigned models with an angular variation of the dust 
density, i.e. 7 = —6.0. The ratio of 7(— 6.0)/7(0.0) drops to 
0.85 for objects with best viewing angles between 40° and 
60° and to 0.5 for objects seen from even higher viewing an- 
gles. The objects with 7 = —6.0 tend to have, on average, 
larger accretion luminosities suggesting that L acc might be 
influencing the dust distribution as well as the torus opening 
angle. 

Lir also increases with increasing Mnust, as larger 
amounts of dust produce stronger IR emission. The mass 
of dust, Moust, increases with redshift as well as with L acc 
(Fig. [6]l . The dependence on the redshift, however, is mainly 
an observational bias: Mnust is computed integrating over all 
dust grains over all volume elements of the torus and there- 
fore depends on Ri n that scales with -y/(Lacc) (from equation 
which increases with increasing redshift (Scott effect). 

In order to make sure all tendencies with redshifts are 
simply the result of the observational bias all flux limited 
samples suffer from, we divided part of the sample to slices 
of Lacc and redshift and computed the mean for MD US t , Lir 
and CF in these bins. The bins were chosen in a way to 
have at least 3 redshift ones per luminosity range. The re- 
sult are given in Table The three quantities show no sta- 
tistically significant trend with redshift in either of the two 
Lacc bins. The average values of the three quantities show 
an increase in lower luminosity bin but can be considered 
constant within the errors. 

The results on the accretion luminosity and the dust 
mass are independent of the choice to impose high optical 
depth models in the fits. The distributions shown in Fig. 
[S]are indistinguishable. This confirms the robustness of the 
derived values of L acc , as further discussed in Section 15.61 



5.3 The torus dimensions 

The inner radius of the torus, i.e. the minimum distance at 
which dust grains can exist, is computed from equation|3]and 
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Table 3. Mj3 ust , Ljr and CF averaged over redshift and luminosity bins. 



bin# 


N 


Lace 






M D ust 


Lir 


CF 






[10 46 erg/sec] 






M Q 


[10 46 erg/sec] 




1 


8 


2.0-5.0 


0.5- 


-1.0 


0.66±1.34 


0.96±0.35 


0.42±0.34 


2 


37 


2.0-5.0 


1.0- 


-1.5 


1.22±1.75 


1.13±0.91 


0.43±0.27 


3 


27 


2.0-5.0 


1.5- 


-2.0 


1.75±2.99 


1.25±0.71 


0.49±0.27 


4 


7 


5.0-10.0 


1.0- 


-1.5 


2.59±3.95 


1.97±1.61 


0.31±0.19 


5 


18 


5.0-10.0 


1.5- 


2.0 


3.88±4.83 


2.90±2.02 


0.52±0.26 


6 


14 


5.0-10.0 


2.0- 


2.5 


2.84±4.50 


2.44±1.24 


0.45±0.28 



for the sample under study typically takes values up to 5 pc, 
with a broad distribution and a peak at ~2 pc. Equation [3j 
however, only takes into account the accretion luminosity as 
a heating source for the dust. And while this is realistic for 
the majority of cases, there are some particular conditions, 
for which the thermal emission from dust itself can become 
an important, non- negligible contributor to its self-heating. 
This this likely to occur in models with steep density pro- 
files (/3 ^ —1) and high optical depths that will give rise 
to so high density at the inner radius that a large amount 
of the dust-emitted radiation is re-absorbed locally and the 
sublimation temperature is reached and exceeded for all the 
dust grains. When this happens during the computation of a 
model, the inner radius is moved to larger values, at steps of 
0.2 pc, and the calculation of the equilibrium temperature 
starts again. This procedure is repeated until the tempera- 
ture at the innermost radius remains below the sublimation 
value. 

As already mentioned, our models use a fixed R ou t /Rin 
ratio, namely of 30 and 100. These cuts are only useful when 
the dust density is assumed constant and does not vary nei- 
ther with the distance to the centre nor with the angle from 
the equatorial plane, i.e. when = 0.0 and 7 = 0.0. When 
= —1.0 though, the dust density rapidly decreases with 
the distance to the centre and the value of R ou t is in practice 
much smaller than that given by the constant R ou t/Rin ra- 
tios. We find the outer radii of the tori to be typically some 
tens of parsecs large but can reach, in a few cases, a few 
hundreds of parsecs. 



5.4 Equatorial optical depth and hydrogen 
column density, Nh 

The hydrogen column density, Nh, is proportio nal to the 
optical depth and is defined as (jFritz et al.l l2006): 

Ndust 

r cq (A) = N H x J2 [Qu(A) + Qfd(A)]af d ^-Nd id (4) 

id=l 

where and Qf d are the absorption and scattering coeffi- 
cients, of the id-th dust grain, respectively, N to t is the total 
number of grain sizes, aid their radius and Nm is the num- 
ber of grain, of a given dimension, calculated with respect 
to hydrogen. The equatorial optical depth, r oq , is in general 
different of that computed along the LoS. 

Fig. [7] shows that high T9.7 models correspond to very 
high equatorial hydrogen column densities, Nh, while allow- 
ing for low rg.7 values of Nh as low as a few times 10 21 are 
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Figure 7. Equatorial hydrogen column density, Nh, distribu- 
tion for all models (red plain line) and high optical depth models 
(black dashed line). 

feasible. In the latter case, the Nh distribution is bimodal 
peaking at the extreme of the values allowed. 

The Nh values as derived from the model may not be 
directly comparable to observed ones (e.g. those computed 
via X-ray observations). In the present study we assume a 
Galactic dust-to-gas ratio, that is likely to vary from one 
galaxy to another. Furthermore, the model value is only ac- 
counting for hydrogen within the torus, while observations 
one takes into account also the gas within the torus inner 
radius, as well as the gas likely present in the outer parts of 
the galaxy as, e.g., diffuse medium. 

5.5 The objects with 70 and 160 um detections 

The objects with 70 um detection represent 25% of the sam- 
ple under study. They have already been included in our 
previous analysis but are worth a closer scrutiny. Their 70 
and 160 um fluxes, as well as their properties, derived from 
the best fit models are shown in Table [4] The distribution of 
their F70/F24 colour versus their 24 pm fluxes, F24, is shown 
in Fig. [8] F70 /F24 usually takes values in the range between 
1 and 8, even for the brightest 24 um sources, depending 
on the relative contribution of the starburst to the IR lumi- 
nosity. For objects with F24 ^ 0.3 mjy there are no 70 pm 
detections because of the 70 pm limiting flux of the survey 
and the 70 pm source definition (see Section [2}. 

Of the 70 objects, 55 had acceptable fits (with \ 2 ^ 16), 
with only one detected at 160 pm. Five of the 55 observed 
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Figure 11. 24 pm flux histogram of the entire sample (open) and 
of the sources with 70 pm counterparts (grayed). 



Figure 8. Observed F70/F24 versus 24 pm flux for the 70 objects 
with 70 pm detections. 



SEDs were reproduced with a single component only, namely 
a torus; the other 50 required a starburst template. Five low 
redshift objects also required a stellar population to fit the 
optical/UV part of their SEDs. As mentioned in Section r3.3l 
two starburst templates were used in the fits, Arp220 and 
M82. Figs. l9l and [TOl shows the observed SEDs (red triangles) 
and the components of the models reproducing them for the 
55 objects with good reduced x 2 and an assigned starburst. 
The torus components are shown in blue dotted lines, the 
starbursts in light green dashed lines and the stellar compo- 
nents in dark green long-dashed lines. The total models are 
illustrated in solid black lines. 

All but one objects were assigned an Arp220-like tem- 
plate to reproduce their 70 pm detection. The reason for 
that is the very bright limiting magnitude of the survey in 
the MIPS Ge bands. An Arp220-like component could make 
the 70 pm bright enough to be observed without contribut- 
ing to the IRAC fluxes while and M82 component would 
also contribute considerably to the IRAC and MIPS 24 pm 
fluxes. 

Another point favouring the conclusion that objects not 
detected at 70 pm could in fact be characterised by a less 
extreme starburst (e.g. M82-like) starburst is the fact that 
even though objects with 70 pm detections are the bright- 
est at 24 pm, not all bright 24 pm sources have a 70 pm 
detection, as seen in Fig. 1111 The opposite, i.e. all bright 24 
pm sources observed at 70 pm, would imply that fainter 24 
pm sources would have different F70/F24 ratios and there- 
fore a different starburst component (M82), but there is no 
observational evidence in favour. 

The contributions of the AGN to the total IR luminos- 
ity, as defined in Section \A. II for the objects with starburst 
components is shown in Fig. [12] (filled circles for the run with 
all T9.7, open circles for the run with high values of T9.7), as 
a function of their 24 pm flux and redshift. The mean dif- 
ference of the AGN contribution derived from the two runs 
is of 3% with a standard deviation of 8%, the results can 
therefore be considered independent of the choice of T9.7S. 
The AGN fraction varies from ~20% up to 100%, with lower 
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Figure 12. The fraction of IR luminosity, as defined in Section 
14.11 attributed to the AGN as a function of their 24pm flux and 
redshift, for the 55 objects with 70 pm detections and \ 2 < 16. 
Filled (open) circles correspond to the run with all (high) 7-9.7. 



redshift objects tending to have higher AGN contribution. 
The increase of the AGN fraction with the 24 pm flux is 
due to the fact that, when the starburst component is an 
Arp220, the 24 pm flux is produced by the torus alone, as 
it corresponds to a very deep absorption feature in the star- 
burst template. This requires a brighter AGN and a lower 
starburst component in order to keep the ratio F70/F24 in 
the observed range. The presence of strong starburst compo- 
nent in many of the cases is an interesting result, when one 
considers the fact that the quasars under study are among 
the brightest in the sky. 

It could be argued that allowing for larger dust distri- 
butions might increase the fraction of the AGN contribution 
to the infrared luminosity or reproduce the 70 micron point 
without the need of a starburst component. In order to ac- 
count for this possibility but also in an effort to reproduce 
the few observed 160 pm datapoints and potentially improve 
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Figure 9. Observed (red triangles) and model (torus: blue dotted line; starburst: dashed light green; stars: long-dashed dark green) 
SEDs for the 55 objects with 70 um detection and minimum \ 2 < 16. 
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Figure 10. Same as in Fig. [9] 
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the fit for the rest of the objects with 70 pm data, we re-run 
the SED fitting ailowing for models with R ou t/Rin = 300, 
a priori excluded from the runs for reasons presented in 
Section T3. II In none of the cases, however, was a torus com- 
ponent alone enough to reproduce the 70 pm point. In fact, 
all objects required a starburst component (again Arp220) 
with a contribution within at most 10% lower than that pre- 
dicted when no large tori were allowed. Furthermore, the 160 
pm points were still not properly reproduced, as all combi- 
nations of torus and starburst templates produced model 
SEDs lying significantly below the observed 160 pm points. 



5.6 Multiple local minimae and degeneracy 

SED fitting in a multi-parameter space usually involves de- 
generacies, i.e. the existence of various combinations of pa- 
rameter values that yield equally good results when repro- 
ducing a set of observed datapoints. Here we are focusing 
on the properties of the AGN emission, and we will hence 
not face this issue with respect to the stellar component 
(both stars and dust heated by star formation). We will ex- 
ploit the fact that the MIR domain is strongly dominated by 
AGN in our sample objects. With respect to this approach, 
the degeneracy problem translates into the fact that there 
are more AGN models, with different parameters, that yield 
equally good -or acceptable in terms of x 2 values- fits, so 
that the properties of the dusty torus can not be unequivo- 
cally assessed. 

To properly deal with this issue, we keep trace of the 
exploration of the model's parameter space, and we analyse 
the 30 best model fits for each object of the sample. For most 
of the objects, even when only the first two solutions are con- 
sidered the discrepancies, measured as Aparam/param(0) 
where param(0) is the value of a given parameter corre- 
sponding to the best fit model, are larger than 50%. In order 
to understand this behaviour one has to take into account 
the influence of each of the model parameters to the global 
model SEDs as well as the fact that not all parameters are 
independent. In fact, most of them are correlated in some 
way. For instance, the LoS and covering factor are tightly re- 
lated: if the covering factor is small, LoS can vary in a larger 
interval; however as the covering factor increases, the LoS 
is restricted in those angles that permit direct view of the 
central source, especially in the case of high optical depths. 

In fact, the parameters that are best constrained are 
L acc and therefore Ri n (from equation [3j and Lir, despite 
the lack of points longward A= 24 pm. Fig. [T3lshows how the 
fraction of objects with standard deviation of the luminos- 
ity over the best fit luminosity, <r(L)/L(0), of ^ 10% (plain 
lines), ^ 20% (dashed lines), and < 50% (dotted lines), 
varies as a function of the first n solutions (n ranges from 
2 to 30). L acc is shown in black while Lir is represented in 
red. 



6 DISCUSSION 

In this work we compare observed and model SEDs of SDSS 
quasars with SWIRE counterparts aiming at constraining 
the model parameters and at quantifying the IR properties 
of bright quasars. The need f or such an appr o ach, t hough 
evident, was also stressed by iGallagher et all ((20071 , who 
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Figure 13. Percentage of objects with cr(L)/L(0) below given 
thresholds, as function of the number of solutions. Black lines: 
Lacc; red lines: Lir. The full, dashed and dotted lines represent 
the following thresholds: ct(L)/L(0) of < 10%, 20%, and 50%, 
respectively. 



studied a sample of 234 SDSS quasars, most of which also 
belong to our sample, trying to quantify the effects of the 
luminosity on the shape of their SEDs in the MIR. Their 
conclusions were solely based on the observed SEDs and 
claimed that comparison with models would constrain phys- 
ical parameters, many of which are dealt with in the present 
study. 

We use a torus model with s moot h dust distribution, 
originally presented in lFritz et al.l l|2006h . even though there 
are conflicting suggestions in the literature as to whether 
tori could be smooth or clumpy. S mooth models have been 
tested on a varie ty of objects (e.g. iGranato fc Danesei 11994 
iFritz et al.ll2006l ) yielding very good results. There is, how- 
ever, ey2dj3nce_J ; lwit_cta be more realistic, 
fe.g. iRisaliti. Elvis fc Nicastrol [20021 ) . nevertheless, the few 
models in the literat ure were tested only on a verage SEDs 
dNenkova et al.ll2002l ) or on individua l objec ts l|Honig et al.l 
I2006T I. Furthermore. iNenkova et"aH (|2002h explicitly ad- 
dressed the issue of suppression of the 9.7 pm silicate fea- 
ture, predicted in emission from smooth models, which was 
never till then supported observationally in type-1 objects, 
but the picture has changed since with Spitzer IRS observa- 



tions of this feature in emission flS icbcnmo rgen et al.l 
Sturm et all 120051 : lHao et alj 120051 : iBuchanan et al 



Ishi et alj|2003 ) 



2005 



2006 



The issue of the influence of the accretion luminos- 
ity on the dust coverage is the source of an ongoing dis- 
cussion in the literature. According to the receding torus 
paradigm (|Lawrencd ll99ll ) the opening angle of the torus 
depends on the power of the central sources, with more pow- 
erful quasars sweeping larger amounts of dust leaving thus 
larger opening angles around them. In this approach, the 
Unified Scheme does not rely uniquely on the orientation 
of the dusty torus but on the dependence of the geometri- 
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Figure 14. Accretion luminosity, L a cc, versus CF: the IllCclll ljacc 
(red squares) decreases in the bins with increasing CF but with 
large dispersions on the values of L acc . 



cal thickness and the optical depths on the central source. 
Evidence for receding tori is often found in studies of radio- 
loud A GN (e.g. Grime s , Raw lings & Willott 2003, 2005). Re- 
cently, iMaiolino et al.l (|2007l ) presented results of a study of 
high luminosity quasars carried out with Spitzer IRS point- 
ing towards decreasing covering factors with increasing lumi- 
nosity. In Section l5.2l we already suggested that the influence 
of L acc on the torus geometry is demonstrated by the larger 
values of L acc generaly occuring in objects with 7 = —6.0. 
In fact, L a cc as derived from the best fit models, also shows 
a slight dependence on CF with the average value of L acc 
in bins of covering factor (red squares) decreases with in- 
creasing covering factor, as seen in Fig. 1141 The dispersion 
on Lacc in each bin, however, is so large that the study is 
inconclusive. 

The study of the sources with 70 pm detections, repre- 
senting 25% of the sample, suggested the presence of a strong 
starburst component in many of the cases. According to the 
best fit models and due to the flux limits of the sample, 
only objects with Arp220-like components were observed in 
70 pm. One could argue that Arp220 is a heavily extin- 
guished starburst, prototype ultraluminous IR galaxy in the 
local Universe and therefore unlikely to be found in high 
redshift bright quasars. In a simpler approach, a black body 
of ~30K could have been used to account for the 70 micron 
detection but we used obs e rved starburst templates instead, 
in the line of iFritz et al.l (|2006l ). However, we did exclude 
from the m odel grid m o dels w ith R ou t/Rm=300, originally 
included in IFritz et al.l l|2006l) . These models imply larger 
dust distributions and lower dust temperatures and could 
have been used to reproduce the 70 pm detections. How- 
ever, they would also correspond to tori with physical sizes 
of several hundred parsecs, sometimes even kpc, where star 
formation is likely to occur. Furthermore, at those distances 
the dust temperature would drop to temperatures typical 



of the dust of starburst or the diffuse dust responsible for 
cirrus emission. 

One of the main points of this work is the use of 
tori models with low equatorial optical depths at 9.7pm, 
T9.7 < 1.0. Our study concludes that the presence of low 
optical depth tori around active nuclei is a possibility, and 
would imply only a minimum of modifications in the current 
picture of the Unified Scheme, namely the possibility of see- 
ing type-1 objects while dust intercepts the line of sight. A 
strong implication of this assumption would be the presence 
of the silicate feature at 9.7 pm in emission even in type-2 
AGN. Such a feature was indeed observed rece ntly in X-ray 
select ed type-2 AGN, using the Spitzer IRS (|Sturm et al.l 
2006). The comparison between results yielded by all T9.7 
models and only high T9.7 ones suggests a flattening of the 
distribution of the covering factor (see Fig.|4j increasing thus 
the apparent estimated ratio between type-2 and type-1 ob- 
jects, as dust enshrouded AGN could still be seen as type-Is 
if seen through a low optical depth medium. In order to re- 
ally probe the Unified Scheme and test the validity of the 
models, on would need to select a complete, volume limited 
sample and study all AGN inside this volume, i.e. of both 
types 1 and 2. The comparison of the tori model parameters 
and the AGN properties between the two types would prove 
or not the validity of the Unification Scheme, since the dis- 
tribution of the various quantities should be the same for 
both types. Furthermore, the SED fitting could provide in- 
formation about the orientation of the sources and finally 
the number of obscured versus unobscured quasars. This is 
the subject of our forecoming work. 

6.1 On the limitations of the models and results 

It has been argued on many occasions that SED fitting tech- 
niques are both powerful and limited. Powerful because pho- 
tometry is and will always be available for samples that are 
orders of magnitudes larger than spectroscopic ones; limited 
because it greatly depends on the quality of the photomet- 
ric data, the model or synthetic SEDs used and the way of 
populating the multi-parameter space, the aliasing and de- 
generacy issues, and the priors one imposes. In the majority 
of cases, the value of SED fitting of large datasets is statis- 
tical, it helps define trends of the population under study 
as a whole, rather than the properties of individual objects, 
unless the sampling of their observed SEDs and the nature 
of the objects themselves allow for it. 

Quasars in general have unmistakable signatures both 
in the UV/optical (in general blue continua attributed to 
the accretion onto the central black hole) and in the IR (IR 
bump due to dust torus emission). Furthermore, the sam- 
ple under study has good wavelength coverage with SDSS, 
2MASS and SWIRE data available, as detailed in Sections 
[2]and[5l helping pin down the wavelength where the falling 
emission from the nucleus meets the rising emission from 
the dust (around 1 micron). 

Since we are working with a pre-constructed grid of 
models, all of the models input parameters take discrete 
values and, in some cases, the model parameter spaces is 
sparsely populated. This is the case, for instance, of the as- 
sumed Rout/Rm ratios, taking two discrete values of 30 or 
100, or the angular variation of the dust density (factor 7 in 
equation [TJ, being either 0.0 or -6.0. Recalculating a model 
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to match each individual object would probably yield bet- 
ter results, however this is computationally very demanding 
and sometimes even impossible, as some of the model pa- 
rameters combination (very high optical depth, T9.7 > 6.0 
and steep density profile s) would lead to a non-convergence 
of the torus model (see iFritz et alj|2006l ). This, of course, 
limits the accuracy of the results and therefore the values of 
the physical parameters derived should be seen as indicative. 
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